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The magnetic susceptibility of the quarter-filled one- 
dimensional extended Hubbard model is calculated using the 
density-matrix renormalization group technique. It is found 
that in the charge gap regime of the model (U > At and 
V > 2t), or in the metallic region with important supercon- 
ductive fluctuations (U < 4t and V > 2t), x(T) displays a 
singularity at T = and an inflection point at low tempera- 
tures that are similar to what occurs in the spin-half quantum 
spin chain. These results, whose accuracy outdoes that of any 
other available technique, are useful data which allow a com- 
parison between theory and experiment in the normal phase 
of the organic conductors. 



During the last three decades, there has been an im- 
pressive collection of experimental results on organic con- 
ductors such as TTF-TCNQ, Bechgaard salts and their 
sulfur analogs. These materials display various low- 
temperature states including superconductivity, itinerant 
and localized spin magnetism, charge and spin density 
wave collective ordering, spin-pcicrls instability, quanti- 
fied nesting and quantum hall effects |lj . It is suspected 
that this rich phenomenology is the consequence of the 
combined influence of low dimensionality and electron 
correlations. These materials are built from large pla- 
nar molecules with 7r-orbitals that are oriented out of 
the molecular plane. Intersite hopping is favored along 
the stacking direction (t a ). The hopping integrals in 
the planar directions are respectively % = t a /\Q and 
t c = ta/500, which lead to strongly anisotropic elec- 
tronic properties. It is generally believed that when 
T > Tx = tb/ir, the physics of these materials is mainly 
one-dimensional. This stems from the fact that there 
is experimental evidence that above Tx, these materi- 
als cannot be described by the Fermi liquid theory. In 
particular, it has been experimentally found that their 
magnetic susceptibility x(T) displays a temperature de- 
pendence which is not compatible with a Fermi liquid. 
Instead of the Pauli behavior of a Fermi liquid, a size- 
able decrease is observed as T is lowered. The Luttinger 
liquid theory (LL) , which plays in one dimension the role 
of the FL, has been proposed to interpret experimental 
results pi . Many of the available LL predictions are valid 
only at very low temperatures, but in this region, the 
3D character reemerges. Thus, one resorts to the direct 
study of model Hamiltonians. 

The ID Hubbard model which belongs to the LL uni- 



versality class is believed to be the appropriate model 
of organic compounds. Although the ground state prop- 
erties of the basic Hubbard model ( V = in equation 
(3) below) were obtained long ago [^0| O], it was only 
recently that its thermodynamic properties were success- 
fully calculated ]12||. There are however experimental 
constraints which imply that the basic Hubbard model 
is not the appropriate model for the organic materi- 
als [||. For instance, the observation of Akp fluctua- 
tions in TTF-TCNQ imply that the correlation func- 
tion exponent K p < 1/2. NMR measurements on the 
(TMTSF) 2 X materials predict that K p = 0.15. These 
values are confirmed by photoemission measurements [|| . 
In the basic Hubbard model, however, K p can not be less 
than 1/2 §. 

The discrepancy between the prediction of the basic 
Hubbard model and the experiments, has led to the sug- 
gestion that dimerization and next nearest neighbor re- 
pulsion V may play a non- negligible role Q [Q. How- 
ever, when these parameters are added, the model is not 
solvable. Knowledge of the strength of correlations is of 
central importance. The correlation parameters can be 
directly extracted from x(^)- Fermi liquid theory and 
band structure calculations predict values that differ by 
an order of magnitude [jj]. Current efforts are oriented 
towards the computation of more reliable values of x(T). 
This quantity is very useful to determine the correla- 
tion parameters. Bourbonnais Q applied a perturbative 
renormalization group (RG) to the g-ology Hamiltonian, 
which for an appropriate choice of parameters is a low 
energy limit of the Hubbard model. He computed x(T) 
and found that its T dependence was important down to 
T = 0. His results are, however, reliable only in the weak 
coupling regime and at low T. Mila and Penc || used 
the world-line quantum Monte Carlo algorithm to com- 
pute x(T) for the dimerized Hubbard model. Aside from 
the usual statistical errors that become important at low 
temperatures, another factor that reduced the quality of 
their calculation was that they had to perform three dif- 
ferent extrapolations on momentum, imaginary time step 
and chain length. A recent analytical RG supplemented 
with a diagrammatic expansion and a QMC calculation 
predicted that x(T) would display a singularity at T = 
and two inflection points at low T for the basic Hubbard 
model ||. These features were first found in the suscep- 
tibility of the spin- half Heisenberg chain ]l3| . These fine 
structures are, however, overlooked by QMC. The rea- 
son is that besides the limited accuracy, the number of 
points computed in the low-temperature region remains 
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small to allow these structures to be resolved. 

In this letter, we show that a thermodynamic density- 
matrix renormalization group (DMRG) algorithm pre- 
sented in an earlier publication |l5f| can lead to an ac- 
curate determination of x(T) for the extended Hubbard 
model. Our results are in excellent agreement with the 
exact results of Ref jl2| for V = 0. We identify two sec- 
tors corresponding to the two regimes of the model. In 
the first, where the Hamiltonian is gapless in both charge 
and spin channels, x(T) increases moderately with V. 
But in the charge gapped regime of the model we observe 
an important rise of x(T). We reach low enough temper- 
atures to clearly see the appearence of an inflection point 
as soon as we enter the regime where the model presents 
a charge gap (V > 2t). At low T, we find that the values 
of x(T) do not smoothly extrapolate to the zero temper- 
ature susceptibility xo- This lead us to conclude that a 
singularity exists at T = 0. 

We apply the thermodynamic DMRG algorithm of Ref 
| p5| Jill . This algorithm is identical to the original T = 
DMRG algorithm . But the computation of the low- 
lying states which was the goal of the T = algorithm 
is now regarded as an intermediary step. In this step, 
we calculate as many target states as necessary to build 
a renormalizcd Hamiltonian which well describes higher 
energy states. For this, we target M states belonging to 
different spin sectors. The reduced density matrix may 
be written as follows, 



(i) 



where pk and ojk are respectively the reduced density 
matrix of the k th state and its weight, p^ is given by the 
wave function relation 
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(i 1 ,i 2 ;i' 1 ,i' 2 ) = i 4 )*j!l(*i> u)- (2) 



As usual, the four indices represent the different blocks 
forming the superblock. The main step of the algorithm 
is the exact diagonalization, by a dense matrix technique, 
of the reduced superblock hamiltonian made of the two 
external blocks. The eigenvectors obtained this way are 
used for the computation of thermodynamic quantities. 
A systematic thermodynamic analysis is made in the fol- 
lowing way. We start by fixing the temperature range. 
x(T) is calculated for initial values of block states m and 
target states M. Then its dependence is analyzed by 
varying m and M. We attribute the same weight to each 
target state, ui = 1/M. The importance of a target space 
in the reduced density matrix is thus given by the num- 
ber of target states fixed for this space. We find that at 
the largest lattice size studied here, a good convergence 
is achieved by attributing 32 states to the S T = sec- 
tor and 4 states to each of the other sectors. We keep 
a total of 300 states in the two external blocks and the 
maximum truncation error is less than p(m, M) < 10~ 4 . 



A detailed analysis of the procedure has been given else- 
where jl6|| . xo = x(T = 0) is calculated by the use 
of the LL relation xo = 2/(NA„(N)), where A a (N) is 
the finite size gap defined as the energy difference be- 
tween the lowest state with St = I and the ground state, 
A a (N) = EiiNj-EoiN). We find that the calculated xo 
is in excellent agreement with Shiba's exact result 1 1 1 . 

The thermodynamic equations of the basic Hubbard 
model were formulated long ago but, it was only re- 
cently that Jiittner et al. succeeded in computing phys- 
ical quantities at finite T. In their study, these authors 
took a different approach which avoided the difficulties of 
the thermodynamic Bethe ansatz equations(BAE). These 
are an infinite set of coupled nonlinear integral equa- 
tions that involve an infinite number of unknown func- 
tions. Actual calculations were possible with only a re- 
stricted basis [0 . Jiittner et al. used the mapping of the 
Hubbard model to the two-dimensional classical Shastry 
model. They then applied the transfer matrix method to 
derive a set of BAE with a finite number of unknowns. 
The resulting nonlinear integral equations are solved nu- 
merically. The DMRG algorithm presented above pro- 
vides a possibility of making an independent verification 
of their results. This algorithm has been quite successful 
in the study of quantum spin chains |l5| and of a mag- 
netic impurity in spin chains |17f| . Its first application 
to fermion systems served for illustrative purposes and 
a detailed analysis of its convergence was made fill . In 
that study, it was mathematically shown that the finite- 
T algorithm is as rigorous as the T = algorithm. The 
only difference is that there are two control parameters, 
m and M , instead of m only. 

For a quantitative analysis with experimental results 
in the normal state of organic conductors, the effects of 
dimerization of the bond-length in the direction perpen- 
dicular to molecular planes should also be considered. In 
this study we consider the extended Hubbard model at 
quarter band-filling, which is the density of many organic 
compounds, 

H = -t^2(cf a c i+la + hc) + U yu^n^ 

i i 

+V} j n i n i+ i. (3) 

i 

The phase diagram of the model(3) has been obtained 
by Mila and coworkers || 0]. They applied the exact 
diagonalization on small clusters as well as analytical 
techniques in various limits to show that the model has 
two regimes. The first regime is a metallic Luttinger liq- 
uid phase that is gapless in both spin and charge chan- 
nels. In the LL phase, the dominant fluctuations are 
spin density waves (SDW); in the strong V region they 
observed important superconductive fluvtuations. The 
second regime is a charge density wave insulator charac- 
terized by an opening of the charge gap while the spin 
degrees of freedom remain gapless. The boundary be- 
tween the two regimes is given by the line going from 
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FIG. 1. DMRG x(T) for U = 1,4,6,8,10,12,14,16 from 
bottom to top. The dots represent the transfer matrix results 
of Ref f|. 

(17, V) = (+oo,2) to (4, +oo). 

Fig. 1 shows that the excellent agreement with BAE 
results seen in the study of quantum spin chains ]lq| can 
also be achieved for the Hubbard model. Our results, for 
U = 1 to 16, were obtained by studying lattices of up 
to N max — 24 and then by extrapolating to N — > +oo. 
They are in excellent agreement with the transfer ma- 
trix calculation of Jiittner et al. The maximum relative 
difference between the DMRG and the transfer matrix 
results is 0.05. This occurs at low temperatures where 
the quality of the DMRG calculations is affected by fi- 
nite size effect. The lowest temperature reached in the 
transfer matrix calculation is T = 0.025. In the DMRG, 
the lowest T, which is roughly given by the finite size gap 
A. a (N max ) depends not only on the size of the system but 
also on the value of U. We find that paradoxically the 
strong coupling regime which is hard to study by analyt- 
ical or QMC methods, is the most favorable since it has 
smaller for a fixed size. We performed an interpo- 
lation between T = and T = A rT (N max ). We discuss 
below the reason of a particular choice of the interpolat- 
ing function. 

A recent study || has predicted that in the weak cou- 
pling regime, (i) the position of the maximum T max of 
x{T) is interaction independent for < U < 4 , (ii) 
the inflection point at T = 0.1 of the non- interacting 
model will survive and will be interaction independent, 
(iii) the appearance of a new inflection point which is 
the consequence of the interaction at lower T, (iv) fi- 
nally, the existence of a singularity at T — in analogy 
with the antiferromagnetic Heisenberg spin-half chain. 
We find that the position of T max is gradually moved 
toward low T as the interaction increases. This DMRG 
result is more likely because in the Hubbard model, the 
weak coupling regime is qualitatively similar to the strong 
coupling regime. Thus one expects the same qualitative 
behavior in the two regions. We can not however directly 
verify the predictions (ii) and (iii) because in the weak 
coupling regime, the variations of x{T) are so small that 
these possible features may be wiped away by extrap- 
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FIG. 2. X (T) for U = 8 and V = 0,1, 2, 3, 4 from bottom 
to top. 

olation of the finite N results. In the strong coupling 
regime, the maximum is shifted to low T, a possible in- 
flection point may appear at even lower T. It is thus 
possible that we have not reached the regime predicted 
by the analytical RG. 

The inclusion of V in the strong coupling regime in- 
duces a clear appearance of the regime described by the 
RG study. In Fig. 2 we show X {T) for [7 = 8 and V = 
toV = 4. When U = 8 and V = 4, we have A a = 0.019. 
The inflection point, which is around 0.025, is clearly 
seen in our calculation. For the interpolation between 
T = and T — A a (N max ), we tried two techniques, 
polynomial and logarithmic. The latter has been tried in 
analogy with the study of the spin-half Heisenberg chain 
|jl3[ . The conformal field theory has shown that at T = 0, 
x(T) displays a logarithmic singularity fi3fl . Since adding 
V leads to an opening of the charge gap, in this regime the 
model bears some analogy with the spin-half Heisenberg 
chain. We find that the fit X {T) = Xo + A/LogT + BT 
is better than a simple polynomial fit. The polynomial 
fit introduces a cusp in x(T). Indeed, we could also fit 
the V — curves with a logarithmic function, but this 
seemed incorrect to us because we had not reached the 
inflection point Tj. The decrease of x(T) f° r T > 1} is 
faster than for T < Tf, it is only in this latter regime that 
a logarithmic fit seems correct. From our calculation, at 
first sight one should relate this low T behavior of x(T) 
to the onset of the insulating phase. But the analysis be- 
low of the weak U regime reveals that these features can 
exist even in a metallic phase. Hence, their occurence in 
the case V = is possible. 

In the weak [7 regime, the model remains metallic for 
all values of V. The charge and spin degrees of freedom 
are both gapless. We observe in Fig. 3 that for [7 = 2, 
x(T) increases with V, as for [7 = 8, although in this 
case the magnitude is smaller. For V = 4, we are in the 
region of the phase diagram where superconducting fluc- 
tuations are important. We find that at low T, the best 
fit to data is logarithmic. Hence, the RG predictions are 
also true in the metallic phase. We conclude by continu- 
ity that this is also valid when V = 0. We were unable to 
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FIG. 3. x(T) for U = 2 and V = 0, 1, 2, 3, 4 from bottom 
to top. 



see these features for V = 0, either because they occur at 
very low T or they are wiped out by the extrapolation. 
The existence of the T — singularity may be the con- 
sequence of the spin-charge separation. It is therefore a 
generic behavior of LL [1^] . 

The magnetic susceptibility is a useful quantity for the 
extraction of correlation parameters U and V from ex- 
perimental data. In Ref || it has been suggested that 
the ratio defined by R = (x( T max) - x( T i)) / x( T max) can 
be used to compare theory and experiment. While we 
will perform an extensive comparison in a longer publi- 
cation, let us here look at the case of (TMTSF) 2 {PF) e . 
In this compound, T max — 275K and Tj = 150 K. A 
rough estimate of R eX p is 0.22. It is also found that the 
resistivity displays insulating behavior. This implies that 
one necessarily has U > 4 and V > 2t. In this regime, 
the best match to the experimental value is obtained for 
U = 8 and V = 4. The corresponding estimate of the 
ratio is Rth = 0.16. This value can be improved either 
by increasing U or V slightly or by including the effect of 
the dimerization in the Hamiltonian. This proliferation 
of parameters renders any precise estimation of U and 
V very cumbersome. It will be necessary to supplement 
the susceptibility comparison with that of another physi- 
cal quantity. The rough estimation made here is in good 
agreement with the values U = 7, V = 2.8 and 0.7 for 
the dimerization found by Mila Q| from the analysis of 
the optical conductivity. A recent reexamination of the 
optical conductivity by Vescoli et al. (2lJ] lead to U = 5 
and V = 2. 

In this study, we have shown that using the finite T 
DMRG algorithm, we can compute accurately the mag- 
netic susceptibility of the extended Hubbard model. We 
found that in the insulating regime and in the region of 
dominant superconductive fluctuations, x{T) displays a 
behavior that is reminiscent of the spin-half Heisenberg 
chain. This is in conformity with the prediction of a re- 
cent RG study. Our results can be used for extensive 
comparisons with experiment. This will be done in a 
longer publication. Our current effort is the calculation 
of the NMR relaxation factor T\. A concurrent analysis 



of x{T) and T\ will provide a better estimation of the 
correlation parameters. 
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